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A sandpile model with stochastic toppling rule is studied. The control parameters and the phase 
diagram are determined through a MF approach, the subcritical and critical regions are analyzed. 
The model is found to have some similarities with directed percolation, but the existence of dif- 
ferent boundary conditions and conservation law leads to a different universality class, where the 
critical state is extended to a line segment due to self-organization. These results are supported 
with numerical simulations in one dimension. The present model constitute a simple model which 
capture the essential difference between ordinary nonequilibrium critical phenomena, like DP, and 
self-organized criticality. 
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I. INTRODUCTION 



The idea of self-organized criticality (SOC) was intro- 
duced to describe the behavior of a class of extended 
dissipative dynamical systems which naturally evolve 
to a critical state, consisting of avalanches propagating 
through the system J]J]. From the very beginning it was 
observed that this new idea has some connections with 
ordinary critical phenomena g. More recently a novel 
mean field (MF) analysis of SOC was presented, which 
pointed out similarities between SOC models and models 
with absorbing states ||. Directed percolation (DP) ||] 
is one of the simplest and most recurrent models with 
absorbing states. Under very general guidelines (local- 
ity, scalar variable, etc.) it has been proposed that a 
wide range of models would fall into DP universality class 
[j§-|§. Alhouh SOC models do not belong to the DP uni- 
versality class they have some connection with DP |) . 

Recently, Tadic and Dhar have shown that a class of 
stochastic sandpile model has some analogy with DP || . 
They studied a directed sandpile model in which unsta- 
ble sites topples with probability p. They observed that 
above a critical threshold p c the system shows SOC, while 
below the system is not critical. The critical probability 
was identified with the threshold for DP in a squared lat- 
tice and the scaling exponents were obtained in terms of 
DP exponents. However, they could not give a detailed 
description of the phase diagram of the model, since their 
analysis was limited to the SOC regime above p c , while 
the state below p c could not be characterized. 

Following the work of Tadic and Dhar we study a class 
of stochastic sandpile models with undirected toppling 
rule. As in their model, sites topples with a probability 
p but now grains are distributed to each nearest neigh- 
bor. In order to provide a theoretical description of the 
model we have generalized the MF theory by Vespignani 
and Zapperi jj| including the control parameter p. In 
this way we obtain the complete phase diagram of the 
model. The existence of a critical probability p c and 



a quasi-stationary state below p c are obtained. Based 
on the MF analysis and on the evolution rules we argue 
that the state below p c is similar to DP, but with different 
boundary conditions. Using this hypothesis we apply the 
scaling theory developed for DP to the present stochastic 
sandpile model. Numerical simulations in one dimension 
support our hypothesis. 

The paper is organized as follows. In section |J we in- 
troduce the dynamical evolution rules for the stochastic 
version of the Bak-Tang-Wiesenfeld (BTW) model. We 
perform a single site MF approximation and determine 
the average densities in the stationary state. It is found 
that the driving rate h and p are the only control param- 
eters and that the system is critical in the line segment 
(h = + , p c < p < 1). Then in section III we argument 
the connection with DP and derive some scaling rela- 
tions. In order to test our predictions we have performed 
numerical simulations in one dimension, the main results 
are presented in section IV. Finally, the summary and 
conclusions are given in section 



II. MF THEORY 

We study a stochastic sandpile model defined as fol- 
lows. An integer variable Zi (height or energy) is as- 
signed to each site of a d-dimensional lattice and energy 
is added to the system at rate h. When a site receives 
a grain and its energy exceeds a threshold z c then, with 
probability p, it relaxes according to the following rules 
Zi — » Zi — g and Zj — ► Zj + 1 at each of g nearest neigh- 
bors. Open boundary conditions are assumed. One may 
call this model non-abelian sandpile model with stochas- 
tic rules. The non-abelian behavior makes it different 
from other stochastic models such as the Manna model 
JlOj l . However, we will simply call it stochastic sandpile 
model. 

The first step towards a comprehensive understanding 
of critical phenomena is provided by mean-field (MF) 
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theory, which gives insight into the fundamental physical 
mechanism of the problem. Thus, we start analyzing the 
stochastic sandpile model through a MF approach. With 
this simple picture we introduce the connection with di- 
rected percolation. 

The first MF theory for sandpile models was intro- 
duced by Tang and Bak [||, and only deterministic top- 
pling rules were considered. Latter Caldarelli et al |ll]] 
generalize this MF theory to sandpile models with certain 
degree of stochasticity in the toppling rules. In particu- 
lar Caldarelli et al studied sandpile models where p is a 
function of Zj, such that < p(z) < 1 below the critical 
threshold z c and p(z) = 1 above. Their MF theory and 
further numerical simulations reveals that the stochastic 
rules, introduced in this way, changes the average of z 
but does not destroy the critical state [jllj . In the sand- 
pile models analyzed by Caldarelli et al p is not a control 
parameter, its average value is determined by the system 
dynamics. In these models the system self-organizes it- 
self to a stationary sate, where (p(z)) is such that the 
system remains in a critical state. On the contrary in 
the stochastic sandpile model considered here p is a con- 
trol parameter. It is expected that for sufficiently small 
values of p the critical state will be destroyed. We have 
therefore to develop a MF theory where p appears explic- 
itly as a control parameter. 

Recently Vespignani and Zapperi ||] have introduced 
a more general framework. As a difference with previ- 
ous theories, their MF approach is not based in some 
particular sandpile model but in general considerations 
which are common to all of them, even extendible to 
other SOC models g|. Within this formalism, SOC ap- 
pears as a special case of nonequilibrium critical phenom- 
ena. They have divided the states each site can assume 
in stable (s), critical (c), and active (a). Stable sites are 
those that cannot become active by addition of energy, 
critical sites are those that become active by addition 
of energy and active sites are relaxing and transfer en- 
ergy to their neighbors. These definitions becomes very 
clear in a deterministic sandpile model. For instance (as- 
suming negligible the probability that a site receives two 
energy grains) sites with z < z c ~ 1 are stable, those with 
z = z c — 1 are critical and those with z > z c are active. 
However, in the stochastic model sites with z > z c top- 
ple with probability p. Again sites with z < z c — 1 are 
stable since they can never become active after receiving 
one energy grain. Nevertheless, those with z > z c — 1 
have not a well defined state. For instance, sites with 
z = z c — 1 may topple after receiving one energy grain 
and, therefore, they cannot be stable. However, they are 
not strictly critical because only a fraction p of them will 
topples after receiving a grain of energy. Hence, in the 
case of stochastic sandpile models the subdivision in sta- 
ble, critical and active does not covert all the possible 
states each site can assume, i.e. p s + p c + p a < 1. 

We divide the states each site can assume in stable 
(s), unstable (it), and active (a) and denote their average 
densities by p s , p u and p a , respectively. The definition of 



stable and actives sites is the same considered by Vespig- 
nani y Zapperi, while unstable sites are now those sites 
that may become active by addition of energy. Under 
these definitions, sites with z < z c — 1 are stable, those 
with z = z c — 1 are unstable and those with z > z c may 
be either unstable or active. Only a fraction p of the 
unstable sites will become active after receiving energy 
and, therefore, are critical sites, i.e. 

Pc=PPu- (1) 

This equality makes the connection between our MF ap- 
proach and that of Vespignani and Zapperi In the 
deterministic limit p = 1 there is no distintion between 
critical and unstable sites. 
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FIG. 1. Schematic representation of the single site approx- 
imation for the stochastic sandpile model. Sites are divided 
in stable (s) unstable (u) and active (s). T mn are the transi- 
tion rates for the state m to n. Transitions which do not take 
place are not represented. 

To study the dynamics of the stochastic sandpile model 
we consider the following Markov process for the average 
densities 

~Qj.P n — ^ ^ ^jnnPm ^ ^ ^nmPn-j (2) 

where T nm are the transition rates from the state n to the 
state m (see fig. |l|). By definition of the model, in one 
step stable sites never becomes active and unstable sites 
never become stable, i.e. T sa = T us — 0. T as — q and 
T au = 1 — q, where q is the fraction of active sites that 
becomes stable after relaxing. In deterministic models q 
may be assumed equal to one || . However in stochastic 
models z may take values large enough compared with z c 
in such a way that an active site main become unstable 
after relaxing, i.e. q < 1. Although an active site with 
z S> z c may remain active due to addition of energy this 
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type of transitions are of second order, they can be ne- 
glected for p a and h small. The transition rates T su and 
T ua depends on the probability per unit time that a site 
receives energy. If p a and h are small then the probabil- 
ity per unit time that a site receives more than one grain 
of energy is negligible, and the probability per unit time 
that a site receive a grain of energy may be approximated 
by 



hi =h + (g - e)p a , 



(3) 



where g — e is the effective number of nearest neighbors 
and e is the dissipation rate per toppling event, an effec- 
tive parameter which account for boundary dissipation. 
If u (p) is the fraction of stable (unstable) sites that be- 
come unstable (active) after receiving a grain of energy 
then T su — uh\ (T ua — phi). Taking into account these 
considerations the system of differential eqs. (|2|) is re- 
duced to 



d_ 

dt Pa 



[1 - (g - e)p c ]p a + p c h + 0{h 2 ,h Pa ), (4) 



d_ 

Of 



Ps = qp a - u(h + gp a )ps + 0(h 2 ,hp a ), 



together with the normalization condition 
Ps + Pu + Pa = 1. 



(5) 



(6) 



Notice that, among unstable sites, only the fraction of 
critical sites p c = pp u contributes to the system dynam- 
ics, the other fraction is only relevant through the nor- 
malization condition in eq. (^) . The system of equations 
is completed by the equation of energy balance 



§- t E=(h- ep a )L d , 



(7) 



where E is the total energy of the system, hL d is the 
average influx of energy and ep a L d the average outflux 
of energy. 



A. Critical state 



In the stationary state (^§jr = 0, ^ =0) from eqs. 



(f§|7|) and (0) we obtain 



Pa = ~, p c =-+0(h), 

e 9 



(8) 



p u = - + 0(h), Ps = PE-J- + 0(h), (9) 

pg pg 



q_ = pg-i 

u p 



+ 0{h). 



(10) 



Comparing this expressions with the ones obtained by 
Vespignani and Zapperi we observe that the average den- 
sities of active and critical sites have the same stationary 
solutions. The differences appear in the density of sta- 
ble and unstable sites, which now depends on the new 
control parameter p. 

For 1/g < p < 1 we have l/<7<p«<l,0<p s < 
(g — 1) / g and, therefore, there is not any inconsistency 
in the stationary solutions obtained above. In this range 
of p, within the MF approach, there is no distinction 
between the critical state of stochastic and deterministic 
sandpile models. The model is critical in the double limit 
h, e — » and h/e — > M and the susceptibility, 



X = 



dp a 

dh 



(11) 



diverges in the critical state. For a small perturbation 
around the subcritical state p a = h/e + Ap a , p c « 1/g 
and from eq. (0) one obtains 



Ap a (t) oc exp(-et/g), 



(12) 



in agreement with the result obtained for deterministic 
models S. 

From eqs. ([ll]) and (|l2| one may think that e is a con- 
trol parameter of the model. However, if the dissipation 
takes place only at the boundary then e will decrease 
with decreasing system size, because the number of ac- 
tives sites in the bulk grows faster than the number of 
active sites at the boundary. Hence, we are just dealing 
with a finite size effect. The statement the system is not 
in a critical state is equivalent to the statement there is 
no critical state in a finite system. If dissipation takes 
place only at the boundary e is not a control parameter, 
it just reflects a finite size effect which at the same time 
is a necessary condition to obtain a stationary state. In 
this sense the criticality here is different from the criti- 
cality at phase transitions where boundary effects always 
disappear in the thermodynamic limit |jj . 

On the contrary, the driving field h is actually a con- 
trol parameter. Since h must satisfy h < e and e — > 
when L — > oo then we must fine tune h to zero in order to 
obtain criticality in the thermodynamic limit. The time 
scale separation becomes a necessary condition for criti- 
cality. Now if we assumes separation of time scales then 
p will be the only control parameter of the model. This 
hypothesis is in general fulfilled in computer simulations, 
where a new grain of energy is added only once there is 
no active site. 



B. Break-down of SOC by stochastic rules 

When < p < 1/g the stationary solutions in eqs. (||- 
p^| ) are no longer valid, because they imply p u > 1 and 
p s < 0. To understand the origin of this inconsistency 
let us analyze the variation of p u and p c with p. In the 
deterministic case p — 1 there is no distinction between 
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unstable and critical sites, i.e. p c = p u — 1/g. How- 
ever, when p < 1 critical sites are a fraction p c = pp u 
of unstable sites. Hence, since p c = 1/g in the station- 
ary state, the system has to self-organize itself increasing 
the average density of unstable sites to p u = 1/pg > 1/g 
and decreasing the average density of stable sites. But 
when p = 1/g we have p u = 1 and p s — and, therefore, 
if we continue decreasing p the system cannot provide 
more unstable sites. Then p c < 1/g and the stationary 
solutions in eqs. ([§f]lO|) break-downs. 

Let us assume that for < p < 1/g the average den- 
sities reache an stationary state, which off course cannot 
be given by eqs. (p[-^0[). From cq. (^) it results that 



Pa 



p c h 



1 - (.9 - e)Pc 



Substituting this expression in eq. (R) one obtains 



d_ 

dt' 



E 



1 



gpc 



i - (g - e)Pc 



-hL° 



(13) 



(14) 



According to eq. (|l]) p c = pp u < p < 1/g, independent 
of p u . Hence, ^ > and the total energy will increase 
linearly with time. Moreover, stable sites are those with 
z < z c — 1 and, therefore, it is expected that after a time 
long enough there will be no stable sites. In this quasi- 
stationary state the energy increases with time and the 
average densities will take the stationary values 



Pa 



hp 



1 - (g - 



P + 0(h), 



p u = l + 0(h), p s = 



q = 0. 



(15) 



(16) 



(17) 



Now it is clear that p is a control parameter of the class 
of stochastic sandpile models analyzed here. 

The susceptibility in this region is given by, assuming 
e<3, 



X 



P 



1 



pg 



(18) 



For a small perturbation around the subcritical state we 
have p a = xh + Ap a and p c w p, then from eq. (||) one 
obtains, again considering e < j, 



Ap a (t) oc exp[-(l -pg)t] 



(19) 



The critical state breaks-down by the stochastic rules, 
once p < p c . To reach the critical state we have to fine 
tune p. Near the critical threshold 1/g the susceptibility 
X ~ (jPc ~ p) and the characteristic time ~ (p c — p)^ 1 
diverges. 

In summary, we found a critical probability p c above 
which the system is in a SOC state, while it is in a sub- 
critical state below. 
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FIG. 2. Evolution of an avalanche in the stochastic sandpile 
model. Empty sites are inactive (stable+unstable) and filled 
sites are active. This figure clearly shows that the avalanche 
evolves in a directed square lattice, where the probability that 
a site is present is p c . 



III. SCALING THEORY 

p c is the probability that a site becomes active after 
receiving a grain of energy. In the absence of an exter- 
nal field it is the probability that a site becomes active 
if one of its neighbors was active in the previous step. 
This problem is equivalent to site directed percolation in 
d + 1 dimensions (d spatial dimensios+time) , p c being 
the probability that a site is present. The only different 
is found in the boundary conditions, while in DP the sys- 
tem is assumed infinite here we deal with a finite system 
with open boundaries. This picture is better represented 
in fig. | 



A. < p < p c 

According to the MF theory below the p c there are no 
stable sites and the fraction of critical sites is given by 
p c = p. The evolution rules are thus that of DP. The 
existence of different boundary conditions may carry as 
a consequence that some scaling exponents result differ- 
ent, however the nature of the phenomena is the same. 
For instance, DP near a wall reveals that the correlation 
length exponents are identical to those obtained in DP in 
an infinite lattice [jl2| , but other exponents take different 
values. This is a consequence of the fact that in DP near 
a wall the avalanches are a subset of the avalanches in DP 
in an infinite lattice. We thus expect a similar behavior 
in the stochastic sandpile model below p c . In this case, 
avalanches starting far from the boundary behaves as in 
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DP in an infinite lattice, while avalanches starting near 
the boundary behaves like the avalanches in DP near a 
wall. Hence, the correlation lengths and the correlation 
length exponents are identical to those of DP in an infi- 
nite lattice, i.e. 



£ii ~ (p c -p)~ 



(20) 



where £j_ and £|| are the spatial and temporal correla- 
tion lengths, respectively, and v± and 1/11 the correlation 
length exponents. 

On the other hand, based on this analogy with DP we 
write the following scaling relation for the average den- 
sity of active sites at site x and time t, given a site was 
active in the origin x — at t = 0, 



p a (x,t)=t T >- 



;f(— - 

Vi 2/Z 'e|| 



(21) 



first introduced by Grassberger and de la Torre in the 
contest of DP [l3| ]. Here r\ is a scaling exponent and z 
the dynamic scaling exponent, as it is usually defined in 
the contest of critical phenomena. Moreover, the prob- 
ability that the avalanche survive up to time t is given 
by 



P(t) 



t- 5 9 



(22) 



where 5 is another scaling exponent. From eq. ( pi|) one 
can derive other scaling laws for the average number of 
active sites n(t), the cluster mass m(t) and the mean 
squared displacement R 2 (t) at time t, resulting 



n(t)=fd d xp a (x,t)=f>f 1 (±), 

R 2 (t) = ^fd d x Pa (x,t)xi = ty*f 3 (±). 

The exponent z is not independent, since t ~ 
eq. ( po|) one obtains 



(23) 



from 



(24) 



and it is therefore identical to that of DP in an infinite 
lattice. Nevertheless, the exponents r\ and S depend on 
the boundary conditions, as it is observed in DP near a 
wall pi. 



B. 



<P< 1 



We assume that the scaling laws in eqs. (|21j-|23|) are 
also valid above p c , but with £|| = In this region 

the dynamical evolution is independent of p and the char- 
acteristic length and time depends only on the lattice size 
L, according to 



However, as it is shown below, for p > p c the global con- 
servation introduces a constraint between the exponents 
rj and z. 

Let us calculate the average flux of energy J(r) outside 
an sphere of radius r, given a grain of energy was added 
at the origin r = at t = 0. The energy flux is pro- 
portional to the gradient of the average density of active 
sites and, therefore, 



d 

J(r) oc r 1 I dt—p a (r,t) 



(26) 



Substituting the scaling relation for p a (r, t) (|2l| ) in this 
expression it results that 

J(r)=r^*- 2 u(j-\ (27) 

Now, conservation implies that J(r) = 1 for r < £j_ ~ L 
and, therefore, 



(1 + 77)2! = 2. 



(28) 



This scaling relation may seem unusual, a more familiar 
expression is obtained if one calculate the mean avalanche 
size 



(s) = / dtn(t) - L {1+ ^ z - L 2 . 
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(25) 



(29) 



This scaling relation was previously obtained by Dhar 
p5| but for a particular sandpile model. We have here 
demonstrated, using scaling arguments, that it holds for 
any sandpile model with global conservation. 



C. Directed models 

In the directed stochastic models there is a preferent 
direction I for the avalanche evolution. This preferent 
spatial direction may be identified with the time direction 
in undirected models and then apply the scaling theory 
developed above. This means that the scaling relations 
in eqs. ( f2lff23] ) are also valid for directed models, but 
now x is a (d — l)-dimcnsional vector in the space of the 
non-preferent directions and t gives the evolution in the 
preferent direction I. 

Another important difference between undirected and 
directed models is the place where topling take places 
during the evolution of the avalanche. In the undirected 
model not only the sites in the avalanche front but also 
sites inside this front may be active, transfering energy 
to their neighbors. On the contrary in directed models 
all active sites are in the avalanche front, i.e. at step t 
all active sites are in the layer I = t. Moreover, actives 
sites in layer I transfer energy only to those neighbors 
in layer I + 1, Hence, in directed models the energy flux 
take places following the preferent direction, while it take 
places in all directions in the undirected case. 
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An inmediate consequence of this difference is that in 
directed models the average outfiux of energy from the 
I = t to the I + 1 layer, which is proportional to the aver- 
age number of active sites in the I — t layer, equals one 
and, therefore, 



n(t) ~ 1. 

Then, form eq. ( p3| ) one obtains that 

T) = 0. 



(30) 



(31) 



The energy balance thus leads to a different constraint 
for the set of scaling exponents (77, z, 8), i.e. the directed 
model is in a diferent universality class. 



IV. NUMERICAL SIMULATIONS AND 
DISCUSSION 

In order to test our predictions we have performed nu- 
merical simulations of the stochastic sandpile model in 
one dimension. We have investigated both regimes of 
the phase diagram, the region similar to DP below p c 
and the SOC region above. In both regions we start with 
a flat pile, i.e. zero height in all sites, and let the system 
evolve to the stationary state. Above p c the stationary 
state is characterized by a constant average energy per 
site, which was taken as the stationary condition. Below 
p c the energy increases linearly with time and, therefore, 
we look for another criterion of stationarity. According to 
the MF theory in the quasi-stationary state below p c we 
have p u = 1, which was taken as the stationary condition. 
In all cases we start measuring after the system reached 
the stationary state. Average where taken over 10 000 
000 avalanches below p c and over 1000 000 avalanches 
above. Below p c we use the lattice size L = 10240, which 
was large enough to avoid finite size effects for the values 
of p considered. Above p c we use p — 0.708 and different 
lattice sizes. 

To obtain an estimate of p c we have calculated the 
correlation lengths in the subcritical state for different 
values of p and fitted the numerical data to the scaling 
laws in eq. (|20|). These magnitudes where computed in 
the simulations using the following expressions 



ill 



Eoo ^ — ^ L 
t = 2-^i = Q P al 

y L tp„ 

V F Pat 



(32) 



where io is the position of the initial active site, t is 
the number of steps measured in the time scale of the 
avalanche and p a i = 1 (p j = 0) in active (unstable) 
sites. 

The log-log plot of the correlation lengths versus p c —p 
is shown in fig. f|. The best fit to the numerical data was 
obtained for 



p c = 0.707 ± 0.002, (0.705485); 
v± = 1.07 ±0.03, (1.0968); 
vn = 1.71 ±0,03, (1.7338). 



(33) 



Enclosed in parenthesis are the series expansion estimates 
for DP in an infinite lattice reported in jl4|] . Within the 
numerical error there is a complete agreement between 
the values reported here and those of DP. 

Then we proceed to determine the exponents 8, rj and 
z from the data collapse plots of P(t), m(t) and R 2 (t), 
using the scaling laws in eq. (^3|). The corresponding 
plots are sown in figs. [||(]. The best data collapse was 
obtained for 



8 = 0.18 ±0.01, (0.15947) 

j] = 0.27 ±0.01, (0.31368) 

z = 1.59 ±0.01, (1.58074) 

u\\ = 1.73 ± 0.01, (1.7338); 
(l + r/)z = 2.02 ±0.02. 



(34) 



From the data collapse we have obtained a better esti- 
mate for v\\ and, using the scaling relation ( |2~i| ) and the 
value of z in (p4), we obtain the better estimate for v±_ 



1.09 ±0.02, (1.0968). 



(35) 



As it was expected the critical probability, the correlation 
length exponents, and z are identical, within the numer- 
ical error, to those reported for DP in an infinite lattice, 
while rj and 8 results different. 

Now let us analyze the numerical simulations in the 
SOC region p c < p < 1. In this case we can obtain an es- 
timate of the critical probability from the divergence of 
the average energy per site (E) near p c . In the SOC 
region {E) reaches a stationary value but it increases 
with time below p c . One thus expect that (E) diverges 
when the system approaches the critical probability from 
above. We observe that the divergence of (E) can be fit- 
ted to the power law dependency (E) ~ (p— p c )~ x , where 
A is a scaling exponent. In fig. we have plot the best fit 
to the numerical data for a lattice size L — 1280, which 
was obtained for 



p c = 0.704 ±0.01, (0.705485), 



(36) 



which is close to the DP value. 

Then we proceed to determine the exponents 8, rj and 
z from the data collapse above p c , using the scaling laws 
in eqs. (^3|) and (^5|). The best data collapse is shown in 
figs. iMwith 



8 = 0.18 ±0.01, (0.15947) 

rj = 0.28 ± 0.01, (0.31368) 

z = 1.57 ±0.01, (1.58074) 
(1 + 77)2 = 2.01 ±0.02 



(37) 



From the comparison of these values with those in eq. 
( |34| ) we conclude that the scaling exponents 8, rj and z 
arc the same above and below p c . Moreover, the scaling 
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relation in eq. (|28|) is in both cases satisfied, although it 
was demonstrating only for p > p c . Hence, there are only 
three independent scaling exponents, v±, v\\ and 6, while 
z and 77 can be determined using the scaling relations 
in eqs. (|24|) and (^). The correlation length exponents 
are identical to those of DP in an infinite lattice while 
S depends on the boundary conditions and, therefore, 
changes the universality class. 
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V. SUMMARY AND CONCLUSIONS 



We have obtained, through a MF analysis, the phase 
diagram of the stochastic sandpile model. There is a 
critical probability p c above which the system is in a 
SOC state, where the correlation lengths diverge in the 
thermodynamic limit. Below p c the system is subcrit- 
ical, it is characterized by a finite susceptibility which 
diverges when the critical state is approached. While the 
stationary state in the SOC state is characterized by a 
well defined average energy per lattice site the subcriti- 
cal state is not completely stationary, since the average 
energy per lattice size increases linearly with time. It 
was then corroborated that the global conservation is a 
necessary condition to obtain SOC in sandpile models. 

Using scaling arguments it was demonstrated that in 
the subcritical region the stochastic sandpile model is 
similar to DP, but with different boundary conditions. 
On the other hand, the scaling theory in the SOC state 
reveals that global conservation introduces a constraint 
among the scaling exponents, generalizing previous re- 
sults obtained for particular sandpile models. We have 
provided a general demonstration of the scaling law (s) ~ 
L 2 . 

Numerical simulations have corroborated the predic- 
tions of the MF and scaling theory. The correlation 
length exponents and the critical probability were found, 
within the numerical error, identical to the estimates for 
DP. However, the existence of different boundary condi- 
tions and conservation law carries as a consequence that 
other exponents result different, changing the universal- 
ity class. 

We must emphasize that the stochastic sandpile model 
is not just another cellular automaton showing SOC, but 
a very nice example to understand the differences and 
similarities between SOC and ordinary noncquilibrium 
critical phenomena. The comparison of the phase dia- 
gram of this model with that of DP reveals the essential 
property of SOC, the insensitive to changes in certain 
"control" parameter, which is off course no more a control 
parameter. While the critical state in DP is restricted to 
a point in the phase diagram, in the stochastic sandpile 
model it is extended through a line segment. 
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FIG. 3. The correlation lengths log£j_ (squares) and log£|| 
(circles) as a function of p in the subcritical state. The lines 
are linear fits to the log-log plot. 
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FIG. 4. Data collapse plot for P(t) in the subcritical stats 
for p = 0.670, 0.688, 0.696 and 0.700. 
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FIG. 5. Data collapse plot for m(t) in the subcritical state 
for p = 0.670, 0.688, 0.696 and 0.700. 




FIG. 6. Data collapse plot for R 2 (t) in the subcritical state 
for p = 0.670, 0.688, 0.696 and 0.700. 
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FIG. 7. Average energy per lattice site as a function of p 
in the SOC state. The continuous line is a linear fit in the 
log-log scale. 
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FIG. 9. Data collapse plot for m(t) in the SOC state for 
L = 160, 320 and 640. 
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FIG. 8. Data collapse plot for P(t) in the SOC state for 
L = 160, 320 and 640. 
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FIG. 10. Data collapse plot for R 2 (t) in the SOC state for 
L = 160, 320 and 640. 
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